Optical coherence elastography to assess biomechanics and detect progression of ocular and other tissues degenerative diseases

ABSTRACT

An excitation force (internal or external) and phase-sensitive optical coherence elastography (OCE) system, used in conjunction with a data analyzing algorithm, is capable of measuring and quantifying biomechanical parameters of tissues in situ and in vivo. The method was approbated and demonstrated on an example of the system that combines a pulsed ultrasound system capable of producing an acoustic radiation force on the crystalline lens surface and a phase-sensitive optical coherence tomography (OCT) system for measuring the lens displacement caused by the acoustic radiation force. The method allows noninvasive and nondestructive quantification of tissue mechanical properties. The noninvasive measurement method also utilizes phase-stabilized swept source optical coherence elastography (PhS-SSOCE) to distinguish between tissue stiffness, such as that attributable to disease, and effects on measured stiffness that result from external factors, such as pressure applied to the tissue. Preferably, the method is used to detect tissue stiffness and to evaluate the presence of its stiffness even if it is affected by other factors such as intraocular pressure (IOP) in the case of cornea, sclera, or the lens. This noninvasive method can evaluate the biomechanical properties of the tissues in vivo for detecting the onset and progression of degenerative or other diseases (such as keratoconus).

This application is a continuation application of and claims priority toU.S. patent application Ser. No. 14/934,663, filed Nov. 6, 2015,entitled “Optical Coherence Elastography to Assess Biomechanics andDetect Progression of Ocular and Other Tissues Degenerative Diseases,”which claims priority to U.S. Provisional Patent Application Ser. No.62/077,561, filed Nov. 10, 2014, entitled “Optical CoherenceElastography to Detect the Onset and Progression of Corneal DegenerativeDiseases,” and U.S. Provisional Patent Application Ser. No. 62/171,043,filed Jun. 4, 2015, entitled “Ultrasound and Optical CoherenceElastography to Assess Biomechanics of Ocular and Other Tissues,” andthe entire contents of these applications are hereby incorporated byreference.

The present invention used in part funds from the National Institute ofHealth (NIH), Nos. 1R01EY022362, 1R01EY014225, 1R01HL120140, U54HG006348and the DOD/NAVSEA, No. PRJ71TN. The United States Government hascertain rights in the invention.

BACKGROUND

This disclosure pertains to a method for assessing changes inbiomechanical properties of ocular and other tissues and for detectingand differentiating tissue stiffness using optical coherenceelastography (OCE).

The changes in viscoelastic properties of the tissues are associatedwith onset and progression of different diseases. Therefore, it isimportant to assess and quantify tissue mechanical properties duringdisease progression and application of different therapeutic procedures.

For example, keratoconus is associated with localized reduced rigidityof the cornea, and the information of the corneal stiffness is useful toprovide improved diagnosis and monitoring of this pathological status.Also, real-time in vivo measurement of the spatial elasticitydistribution with microscopic scale in the cornea could lead to adaptivemechanical modeling of the individual corneal structure which isextremely important to prevent over-corrections, under-corrections andectasia from refractive surgeries, such as LASIK, and to furtheroptimize the laser ablation procedures.

Structurally degenerative diseases such as keratoconus can significantlyalter the stiffness of the cornea, directly affecting the quality ofvision. Keratoconus can pathologically decrease the stiffness of thecornea, leading to a loss in the quality of vision. Detecting changes inthe biomechanical properties of ocular tissues, such as stiffness of thecornea, can aid in the diagnosis of these structurally degenerativediseases.

As another example, changes in mechanical properties of crystalline lensplay an important role in the development of presbyopia, which is theprogressive, age-related loss of accommodation of the eye. Results ofex-vivo studies have shown that the stiffness of crystalline lensesincreases with age for animals and humans. The increase in lensstiffness is generally believed to be responsible for the progressiveloss of the ability of the lens to change shape, leading to presbyopia.However, current understanding of the mechanical properties of the lens,its changes with age, and its role in presbyopia is limited, in part dueto a lack of technology that allows measurements of the mechanicalproperties of the lens in situ and in vivo. The location of thecrystalline lens inside the eye makes it challenging to measure itsmechanical properties in vivo or in situ (i.e., inside the globe).

UV-induced collagen cross-linking (CXL) is an emerging treatment thateffectively increases corneal stiffness and is applied clinically totreat keratoconus. The effectiveness of this treatment may be analyzedby measuring the corneal stiffness both before and after treatment.However, measured corneal stiffness is also influenced by intraocularpressure (IOP). Therefore, experimentally measured changes in cornealstiffness may be attributable to the effects of CXL, changes in IOP, orboth. There is a possibility that a cornea, particularly after treatmentwith CXL, may still be structurally weakened by keratoconus, yet have a“normal” measured stiffness due to an elevated IOP. Current techniquesare not able to measure the true IOP in vivo without consideration ofthe effect of the biomechanical properties of the cornea. Distinguishingcorneas that have the same measured stiffness but are at different IOPsis still a challenge.

Elastography is an emerging technique that can map the local mechanicalproperties of tissues. Ultrasound elastography (USE) and magneticresonance elastography (MRE) have experienced rapid development duringthe past few years as diagnostic tools. One common principle of thesetechniques is correlating tissue deformation caused by the externalmechanical excitation to tissue elasticity. In previous studies,acoustic radiation force was applied to a microbubble created bylaser-induced optical breakdown in the lens. The displacement of themicrobubble was measured by ultrasound and used to evaluate lenselasticity. However, such approach is required to produce microbubbleswithin the lens. The basic feasibility of using Brillouin microscopy tomeasure the lens bulk modulus both in vitro and in vivo has beenexplored. Brillouin microscopy can be implemented using simpleinstrumentation, but it has a relatively slow acquisition time. There isalso uncertainty on how to correlate Brillouin shift (modulus) to theclassical mechanical description of the tissues (e.g. Young modulus).USE and MRE can assess mechanical properties of tissue but therelatively low spatial resolution of USE and MRE is still a criticallimitation for certain applications, particularly for ocular tissues andalso for measurements at the cellular level.

What is needed, therefore, is an improved, non-invasive and highlysensitive method to assess the mechanical properties of the ocular andother tissues with high resolution and sensitivity.

SUMMARY

The present disclosure relates generally to methods and systems forassessing the biomechanical properties of tissues non-invasively, to amethod using optical coherence elastography (OCE) for detecting tissuestiffness, such as corneal or lens stiffness, and for differentiatingsamples having similar measured stiffness as a result of otherinfluences, such as intraocular pressure (IOP). The methods describedherein for tissue biomechanical quantification are demonstrated in thecase of the cornea and crystalline lens but generally applicable for allsoft and hard tissues in the body.

Optical coherence elastography (OCE) is capable of direct and highresolution assessment of mechanical properties of tissue and, therefore,overcomes the limitations of previously-used techniques. OCE employshigh-resolution optical coherence tomography (OCT) to detect the sampledeformation induced by an external force. In comparison to ultrasoundelastography (USE) and magnetic resonance elastography (MRE), OCE isable to provide superior spatial imaging resolution, faster acquisitionspeed, and greater displacement sensitivity.

In one aspect, this disclosure relates to a method for quantifyingbiomechanical properties of a tissue, comprising: producing an externalor internal force to stimulate localized deformation on a surface of thetissue; using an optical coherence tomography (OCT) or any otherlow-coherence interferometry subsystem to measure an induceddisplacement profile resulting from the localized deformation on thesurface of the tissue; and quantifying the biomechanical properties ofthe tissue based on the analysis of the induced elastic wave (i.e.temporal or spatial displacement profile analysis) using an algorithm.

The present system combines a pulsed ultrasound system (or any otherexcitation methods such as air-puff, laser pulse, etc) capable ofproducing a force on the tissue surface and a phase-sensitive OCT systemfor measuring the tissue displacement caused by the force. The systemallows for a non-invasive and highly sensitive method to assess themechanical properties of the tissue in vivo.

The present noninvasive measurement method also uses phase-stabilizedswept source optical coherence elastography (PhS-SSOCE) to distinguishbetween corneal stiffness attributable to disease or UV-induced collagencross-linking (CXL) and IOP effects on measured corneal stiffness.Optical coherence elastography (OCE) is an emerging noninvasivetechnique used in recent years that can map the local biomechanicalproperties of tissue. Similar to ultrasound elastography (USE) andmagnetic resonance elastography (MRE), OCE is usually comprised of anexternal loading component that produces displacements within thetissue. In OCE, imaging this tissue displacement is performed withoptical coherence tomography (OCT), which has superior spatialresolution compared to USE and MRE. From the velocity of an inducedelastic wave (EW), or stress-strain curve measured by OCE, tissueelasticity can be quantitatively estimated.

The present method compares the displacement amplitude attenuation,elastic wave speed, dispersion of the elastic waves, and naturalfrequency of the vibrations of a focused air-pulse induced elastic wavein corneal tissue. The damping speed of the displacement amplitudes ateach measurement position along the wave propagation are compared forthe different materials. This noninvasive method has the potential todetect the early stages of ocular diseases such as keratoconus or to beapplied during CLX procedures by factoring in the effects of IOP onmeasured corneal stiffness.

FIG. 1 shows a flow chart illustrating an example of a methodencompassed by the present disclosure. Generally, in the first step ofthe method, an elastic wave is generated in two corneal tissue samples,preferably by a focused air pulse. The elastic wave is detected usingPhS-SSOCE at different positions along the wave propagation. If themeasured elasticity is different, then the two samples can bedifferentiated without additional measurement. If the measuredelasticity is the same, or at least so close that additional measurementis warranted, then additional steps are taken. The wave displacement isnormalized by dividing the amplitude at each measured position by theamplitude at the excitation position. Then, the attenuation of thenormalized elastic wave displacement is compared. The sample that hasthe faster attenuation has a higher viscosity than the other sample,allowing for the samples to be differentiated.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 shows a flow chart describing an exemplary method fordifferentiating samples based on elastic wave measurements andphase-stabilized swept source optical coherence elastography(PhS-SSOCE).

FIG. 2 shows an example of (a) a schematic of the co-aligned ultrasoundand optical coherence elastography (US-OCE) system and (b) a typicaloptical coherence tomography (OCT) image of a crystalline lens.

FIG. 3 shows (a) elastic wave (EW) velocities of 14.0% gelatin (w/w) and1.1% (w/w) agar phantom samples measured by PhS-SSOCE; and (b) acomparison of estimated and measured Young's moduli of the 14.0% gelatinand 1.1% agar phantoms.

FIG. 4 shows (a) a comparison between the normalized elastic wavedisplacement amplitude attenuation curves of the 14.0% agar phantom attwo different excitation pressures; and (b) a comparison of thenormalized elastic wave displacement amplitude attenuation curves of the14.0% gelatin and 1.1% agar phantoms.

FIG. 5 shows a comparison of the normalized elastic wave displacementamplitude attenuation curves and exponential fit of an untreated (UT)and UV-induced collagen cross-linked (CXL) porcine cornea at the samemeasured stiffness but different intraocular pressures (IOPs).

FIG. 6 shows (a) a schematic and mathematical description of a rabbitlens shape and (b) a transformed cylinder with the same volume as therabbit lens.

FIG. 7 shows (a) temporal displacement profiles of young (n=3) andmature (n=4) lenses and (b) maximum displacements of the young (n=3) andmature (n=4) lenses.

FIG. 8 shows (a) a recovery process fitted by y(t)=A(1+bt)e^(−ωt) and(b) the natural frequencies ω for the young (n=3) and mature (n=4)lenses.

FIG. 9 shows (a) typical experimentally measured and (b) theoreticallycalculated displacements on the lens surface for young and maturelenses.

FIG. 10 shows (a) Young's modulus and (b) shear viscosity modulus ofyoung (n=3) and mature (n=4) lenses estimated based on a model of theviscoelastic layer.

FIG. 11 shows (a) uniaxial mechanical tests and fitted stress-straincurves for typical young and mature rabbit lenses and (b) thedistribution of the Young's modulus E_(0.1) at strain ε=0.1 for theyoung (n=9) and mature (n=4) rabbit lenses.

DETAILED DESCRIPTION OF PREFERRED EMBODIMENTS

The present disclosure relates to methods utilizing optical coherenceelastography (OCE) to detect tissue stiffness and to distinguish theeffects of other factors that might affect tissue stiffness. Previousstudies have demonstrated that OCE is feasible for quantitativelyassessing the elasticity of a sample. Preferably, the method is used todetect corneal stiffness and to evaluate the presence of cornealstiffness even if it is affected by intraocular pressure (IOP). Thepresent method uses phase-stabilized swept source optical coherenceelastography (PhS-SSOCE) and can distinguish untreated (UT) andUV-induced collagen cross-linked (CXL) corneas of the same measuredstiffness but at different IOPs. This noninvasive method can evaluatethe biomechanical properties of the cornea in vivo for detecting theonset and progression of corneal degenerative diseases such askeratoconus.

Generally, the present method is for measuring tissue biomechanicalproperties (e.g. stiffness) and for differentiating tissue samples usingoptical coherence elastography. In an exemplary embodiment, a first stepis inducing elastic waves in the tissue samples, followed by detectingproperties of the waves using interferometry (low-coherenceinterferometry or optical coherence tomography is preferred) atdifferent measurement positions along the waves (or at the same positionfor temporal analysis of the elastic wave displacement profile), whereinthe detected properties include measured wave velocities and measuredwave displacement amplitudes. A next step is determining elasticities ofthe tissue samples using the measured wave velocities and thendifferentiating the tissue samples having different measured wavevelocities. For those samples having similar measured wave velocitiesand needing further differentiation, a next step is normalizing themeasured wave displacement amplitudes to produce normalized wavedisplacement data. The normalized wave displacement data is then used toidentify tissue samples having faster wave attenuation and slower waveattenuation. The tissue samples having faster wave attenuation are thenclassified as tissue sample having increased viscosity and reducedstiffness and the tissue samples having slower wave attenuation areclassified as tissue samples having reduced viscosity and increasedstiffness.

In preferred embodiments, the tissue samples are ocular tissue samples.In addition, the step of inducing elastic waves may be by directing afocused air-pulse on the tissue samples. The step of determiningelasticities of the tissue samples is preferably by calculating Young'smodulus using the measured wave velocities. The step of normalizing themeasured wave displacement amplitudes is preferably by dividing themeasured wave displacement amplitudes at the different measurementpositions along the waves by the measured wave displacement amplitude atan excitation position. The step of using the normalized wavedisplacement data to identify tissue samples having faster waveattenuation and slower wave attenuation preferably involves calculatinga customized ratio between normalized data collected at variouspositions for two samples, as described more fully below. This step canbe repeated for additional pairs of samples.

The present disclosure also relates to a co-aligned focused ultrasound(or any other excitation methods such as air-puff, laser pulse, etc.)and phase-sensitive optical coherence elastography (US-OCE) system forassessment of biomechanical properties of tissues, including in situ andin vivo.

An example demonstrated herein is the use of ultrasound radiation as theexcitation force and crystalline lens as the target. However, any typeof the excitation force (external or internal) or any ocular or othertissues could be used with this method.

The exemplary US-OCE system includes an ultrasound delivery sub-systemthat emits acoustic force to stimulate localized deformation on thecrystalline lens surface and an optical coherence tomography (OCT)sub-system to measure the induced displacement profile. The US-OCEsystem is preferably used in conjunction with a data analyzingalgorithm. Based on measurements from the US-OCE system, the dataanalyzing algorithm can quantify the mechanical parameters of thecrystalline lens. Any suitable computer or data processor programmedwith the data analyzing algorithm can be used to make thesecalculations. The mechanical parameters include parameters such asdisplacement amplitude, natural frequency, Young's modulus, and shearviscosity of the crystalline lens. The kinematical differential equationcan be used to describe the lens's relaxation process starting from themaximum displacement point:

$\begin{matrix}{{{\frac{d^{2}{y(t)}}{{dt}^{2}} + {2{\xi\omega}\frac{{dy}(t)}{dt}} + {\omega^{2}{y(t)}}} = 0},} & (1)\end{matrix}$where ξ=c/(2√{square root over (mk)}) is the damping ratio andω=√{square root over (k/m)} is the undamped natural frequency of thedynamic system. The analytical solution of equation (1) isy(t)=A(1+bt)e ^(−ωt)  (2)when ξ=1; equation (2) can be used to calculate natural frequency. Ananalytical solution of the spectral component of the verticaldisplacement in the frequency domain can be derive:

$\begin{matrix}{{{Y_{z}\left( {r,z} \right)} = {- {\int\limits_{0}^{\infty}{\alpha^{2}{{J_{0}\left( {\alpha\; r} \right)}\left\lbrack {{A_{1}e^{{- \alpha}\; z}} - {A_{2}e^{\alpha\; z}} + {\alpha\left( {{B_{1}e^{- {\partial z}}} + {B_{2}e^{\partial z}}} \right)}} \right\rbrack}d\;\alpha}}}},,} & (3)\end{matrix}$where J₀ and J₁ are Bessel functions of order 0 and 1, respectively.Using the analytical solution of the forward problem (3), reconstructionof Young's modulus and shear viscosity was posed as a minimizationproblem

FIG. 2(a) shows a schematic of an example of the co-aligned US-OCEsystem. The setup consists of the Ultrasound (US) excitation system andOCT detection system. The function generator and the amplifier are usedto drive the ultrasound transducer. FIG. 2(b) shows an example of atypical OCT image of a lens.

For experimental validation, the biomechanical properties of rabbitcrystalline lenses were assessed in situ by using a US-OCE system.Experiments were performed on the lenses of young and mature rabbits insitu (lens located inside an eye globe). Both the maximum displacementand the relaxation rate of the displacement were analyzed. Also, amodel-based reconstruction was applied to quantify the viscoelasticproperties of the lenses. For validation, uniaxial mechanicalcompression tests (considered as a “gold standard”) were conducted onthe same young and mature rabbit lenses. The US-OCE system, whichcombines acoustic radiation force excitation and phase-sensitive OCT,was demonstrated as a promising tool for noninvasive assessment of thechanges in the biomechanical properties of the crystalline lens in situ.The high displacement sensitivity of phase-resolved OCT detectionenables the measurement of sub-micron displacements on the lens surface,which is critical for in vivo study as it allows for the application ofa minimal acoustic radiation force to induce a detectable displacementand minimizes the potential ultrasound damage to the eye. In addition,the high spatial resolution of OCT allows highly-localized investigationof the mechanical properties of the lens.

Example 1. Measurement and Validation

A PhS-SSOCE system was utilized which consisted of a focused air-pulsedelivery system and a phase-stabilized swept-source OCT (PhS-SSOCT)system. A short duration focused air-pulse was expelled through anelectronic solenoid controlled air gate and induced an elastic wave inthe sample. A pressure gauge provided air source pressure control andmeasurement. The localized air-pulse excitation was positioned with a3-D linear micrometer stage. The PhS-SSOCT system was comprised of abroadband swept laser source (HSL2000, Santec, Inc., Torrance, Calif.)with central wavelength of 1310 nm, bandwidth of ˜150 nm, scan rate of30 kHz, and output power of ˜29 mW. A-scan acquisition was triggered bya fiber Bragg grating. The axial resolution of the system was ˜11 μm inair. The experimentally measured phase stability of the system was ˜16mrad, which corresponded to ˜3.3 nm displacement in air. Bysynchronizing the focused air-pulse with consecutive M-mode images, theelastic wave velocity and a two dimensional depth resolved elasticitywere calculated.

A validation study was initially conducted on 14.0% gelatin (w/w) and1.1% agar (w/w) phantom samples (n=5 for each type) with the samecylindrical dimensions of diameter D=33 mm and height H=11 mm. As shownin FIG. 3(a), the EW velocity, c, measured by PhS-SSOCE in the gelatinsamples was 3.76±0.2 m/s, which was very similar to the EW velocity inthe agar samples of 3.64±0.3 m/s. The acoustic surface wave formula(equation 1, below) was used to estimate the Young's moduli of thesamples, where the density, ρ=1.02 kg/m³ and Possion ratio, υ=0.49.

$\begin{matrix}{E = \frac{2{\rho\left( {1 + \upsilon} \right)}^{3}c^{2}}{\left( {0.87 + {1.12\upsilon}} \right)^{2}}} & (1)\end{matrix}$

As shown in FIG. 3(b), the Young's modulus for the 14.0% gelatin and1.1% agar phantoms obtained by the analytical model was 48.7±9.2 kPa and46.6±8.2 kPa, respectively. Uniaxial mechanical compression tests (Model5943, Instron Corp., MA) were conducted on the phantoms for elasticityvalidation. The measured Young's modulus was 47.6±5.3 kPa for the 14%gelatin and 44.9±6.6 kPa for the 1.1% agar phantoms as shown in FIG.3(b). These results demonstrated that the 14.0% gelatin sample and 1.1%agar phantoms were of similar stiffness confirmed by both analyticalmodel and uniaxial compression tests.

To compare the damping characteristics between any two normalizeddisplacement amplitude attenuation curves of the elastic waves, acustomized ratio, r, was used, where ND1i and ND2i were the normalizeddisplacement of the induced elastic wave at the i^(th) measurementposition for samples 1 and 2, respectively.

$\begin{matrix}{{r_{({{ND}_{1}\text{/}{ND}_{2}})} = {{{{mean}\left( r_{i} \right)} \pm {{{std}\left( r_{i} \right)}\mspace{14mu}{with}\mspace{14mu} r_{i}}} = \frac{{ND}_{1i}}{{ND}_{2i}}}},} & (2)\end{matrix}$Displacement amplitudes were normalized by dividing the elastic wavedisplacement amplitude at each measurement position by the displacementamplitude at the excitation position. If r was significantly greaterthan 1, the displacement in sample 2 damped faster than in sample 1. Ifr was significantly less than 1, sample 1 damped faster than sample 2.If r was close to 1, the damping was similar in both samples.

This ratio was first calculated for the same 14.0% gelatin phantom toexamine the effects of different initial position displacements bychanging the focused air-pulse pressure on the sample to 11 and 22 Pa.The normalized displacement attenuation curves are shown in FIG. 4(a)with the ratio r_((22/11))=0.95±0.12, which was very close to 1. Asanticipated, this indicated that the initial displacement amplitude didnot affect the damping speed of the elastic wave.

This ratio was then calculated to compare the gelatin and agar phantoms.As shown in FIG. 4(b), the normalized displacement in the agar phantomswas higher than in the gelatin phantoms at the same scan position. Byusing formula (2), r_((agar/gelatin))=1.56±0.47, which demonstrated thatthe 14% gelatin damped faster than the 1.1% agar. This result was inagreement with previous findings that gelatin has higher viscosity thanagar, which corresponds to faster damping. Therefore, these comparisonsshowed that this method could be successfully utilized to distinguishtwo materials of similar stiffness.

Example 2. Corneal Stiffness

To induce a similar measured corneal stiffness in untreated (UT) andUV-induced collagen cross-linked (CXL) porcine corneas, the IOP of thewhole eye was controlled by a custom-built controller comprising of apressure transducer and micro-infusion pump connected in a feedbackloop. The elastic wave was measured in a porcine cornea by the PhS-SSOCEsystem before and after UVA-Riboflavin induced CXL. Elastic wave (EW)measurements were taken at IOPs from 15-35 mmHg with 5 mmHg increments.The EW velocities of the elastic wave in the UT and CXL corneas at thevarious IOPs are presented in Table 1 below. It can be observed thatbefore CXL, the EW velocity of the cornea at IOP=30 mmHg was calculatedas c=3.6±0.4 m/s. After CXL, the EW velocity was 3.6±0.1 m/s at IOP=20mmHg. Therefore, based on the EW velocity, it might appear that thestiffness of the cornea is the same in those two occurrences.

TABLE 1 UT CXL IOP (mmHg) EW velocity (m/s) EW velocity (m/s) 15 1.5 ±0.1 2.7 ± 0.1 20 2.3 ± 0.1 3.6 ± 0.1 25 3.0 ± 0.3 4.0 ± 0.1 30 3.6 ± 0.44.2 ± 0.2 35 3.7 ± 0.4 4.7 ± 0.5

After normalizing the elastic wave displacement amplitudes, the dampingfeatures of the elastic wave over the measurement positions wereanalyzed (FIG. 5). Based on formula (2), the ratio of r_((CXL/UT)) wascalculated as r=2.61±0.95, indicating that the damping speed of thecornea had significantly decreased after the CXL treatment. In addition,the normalized displacement attenuation curves were fitted by y=ae^(bx)in which the parameter b was treated as the damping speed. According tothe fitted results, the damping speed of the UT cornea (b_(UT)=−0.031mm⁻¹) was almost twice the damping speed of the CXL cornea(b_(CXL)=−0.017 mm⁻¹), which confirmed that the damping speed decreasedafter CXL treatment. One possible reason for this result is that the CXLtreatment is a procedure which displaces water from the cornea tissue.The UT cornea contains more water, which is responsible for a higherviscosity. Therefore, the elastic wave damps faster in the UT corneathan the CXL cornea.

Previous studies have shown that viscosity is negatively correlated withmeasured corneal stiffness, indicating that the CXL cornea has lowerviscosity than the normal one, which corroborates with these results.

Example 3. Experimental System and Sample Preparation

An example of a co-focused and co-localized ultrasound and opticalcoherence elastography system, termed US-OCE, was developed by combiningultrasound excitation with spectral-domain OCT, as schematically shownin FIG. 2(a). A single element ultrasound transducer (model ISO305HP,CTS Valpey Corporation, MA) with a diameter of about 12.7 mm and a focallength of about 19 mm operating at 3.7 MHz center frequency was employedin the system. A sinusoidal wave with a frequency of 3.7 MHz generatedby a 50 MHz function generator (model 80, Wavetek Ltd., CA) was gated bya pulse with duration of 1.1 ms to produce a one-time burst. The drivingsignal for the ultrasound transducer was amplified using a 40 dB poweramplifier (model 350L, Electronics & Innovation Ltd., NY). The acousticradiation force from the ultrasound wave was used to remotely perturbthe anterior surface of the crystalline lens through the cornea and theaqueous humor of the eye.

In the phase-sensitive OCT system, a superluminescent laser diode (modelS480-B-I-20, Superlum Diodes Ltd., Ireland) was utilized as the lightsource with a central wavelength of about 840 nm and a bandwidth ofabout 49 nm. The laser beam was separated and directed to the referenceand the sample arms of a Michelson interferometer. The interference ofthe combined light from these two arms was detected using ahigh-resolution spectrometer comprised of a grating and a line-scanningCCD camera (model L104-2k, Basler, Inc., Germany). The A-lineacquisition rate of this system was set to 25 kHz during theexperiments. A full width at half maximum (FWHM) of the transverseGaussian profile of the OCT beam at the imaging focal plane was about 8μm. The system's phase stability measured by the interference signalfrom the reflection of the two surfaces of a glass slide in the samplearm was ˜4 milliradians. However, any OCT or other low-coherenceinterferometry system, which can measure nanometer to micrometeramplitude displacements, can be used with this method.

A custom-built transducer holder was used to securely attach theultrasound transducer to the OCT objective lens. The co-alignment of thefocal zone of the ultrasound beam and the OCT imaging beam was achievedby aligning the mounted ultrasound transducer with a needle tip.Acoustic radiation force excitation and OCT M-mode imaging (rapidlyrepeated A-scans at the same location) were synchronized by acomputer-generated triggering signal.

Eyes from three young (2-3 months old) and four mature (over 6 monthsold) rabbits (Pel-Freez Biologicals, LLC, AR) were used in this study.Immediately after enucleating, the globes were placed in a 1×phosphate-buffered solution (PBS, Sigma-Aldrich Inc., MO) and shippedovernight over a dry ice (without freezing). All experiments wereperformed immediately after receiving the eyes. During the experiments,the entire eye globes were kept in the saline at room temperature tominimize any change in the tissue properties. The sample was positionedin a custom-designed eye holder to prevent motion during the experiment.

The surface of the crystalline lens was placed at the co-aligned focalzone of the US-OCE system. The axis of the OCT beam and, therefore, thedirection of the measured displacements was orthogonal to the lenssurface. However, the ultrasound transducer was placed at an angle ofabout 45° relative to the OCT sample beam, so both axial (i.e., verticalor along the axis of the OCT beam) and transverse components of acousticradiation force were generated. Excitation with the acoustic radiationforce produced a perturbation on the lens surface, resulting in adisplacement of the lens surface. The displacement of the apex of thecrystalline lens as shown in FIG. 2(b) was measured by thephase-sensitive OCT system. Only the axial displacement at the apicalposition of the lens surface was measured by the US-OCE system. Ingeneral, the measured axial displacements are related to the stiffnessgradient inside the lens, including the cortex, nucleus, and the lenscapsule. The measurements only reflect the elastic properties in thearea of the ultrasound excitation where higher elasticity leads tosmaller maximum displacement. Since the ultrasound pulses were deliveredthrough a part of the cornea, outside of the optical path of the OCTbeam, displacement of the cornea in the OCT image was minimal and didnot contribute significantly to the measured signal from the lenssurface.

During the experiment, the distance between the sample and theultrasound transducer was held constant. Therefore, the acousticradiation force applied on the lens surface can be considered asapproximately the same for all the samples, which eliminates theinfluence from the magnitude of the acoustic radiation force on theamplitude of the displacement on the lens surface.

Example 4. Modeling

Kinematic Model of the Relaxation Process

Under an external acoustic radiation force, the movement of the tissuein the focal zone of US-OCE, shown as a dot in FIG. 2(b), can beseparated into forced motion in the response to the acoustic radiationpush, and relaxation motion that occurs after the acoustic radiationforce is removed and when the external forces are zero.

The following simplified kinematical differential equation can be usedto describe the lens's relaxation process starting from the maximumdisplacement point:

${{{m\frac{d^{2}{y(t)}}{{dt}^{2}}} + {c\frac{{dy}(t)}{dt}} + {{ky}(t)}} = 0},$where m is the equivalent mass, c is the viscosity coefficient and k isthe equivalent spring stiffness. To understand the basic characteristicsof the equation, two parameters, ξ and ω, are introduced whereξ=c/(2√{square root over (mk)}) is the damping ratio and ω=√{square rootover (k/m)} is the natural frequency of the dynamic system. The equationthen becomes

${\frac{d^{2}{y(t)}}{{dt}^{2}} + {2{\xi\omega}\frac{{dy}(t)}{dt}} + {\omega^{2}{y(t)}}} = 0.$

The analytical solution of the second equation is related to the valueof ξ as:y(t)=A(1+bt)e ^(−ωt) when ξ=1; and  (a)y(t)=e ^(−ξωt)(A cos ω_(D) t+B sin ω_(D) t) with ω_(D)=ω√{square rootover (1−ξ²)} when 0<ξ<1.  (b)Here A, b, and B are the parameters to be determined. According to theexponent forms of the solution of the second equation, ω can also bedescribed as the relaxation rate, which corresponds to the rate of theexponential-type displacement recovery process.

Model for a Viscoelastic Layer

To quantitatively evaluate age-related changes in the viscoelasticproperties of the rabbit lenses, a model-based reconstructive approachbased on the deformation of a homogeneous viscoelastic layer in responseto an acoustic radiation force of short duration was considered. In thisapproach, tissue is modeled as an incompressible viscoelastic (Voigtbody) layer. An acoustic impulse is considered as an axisymmetric forceapplied to the upper surface of the medium in the direction of thez-axis of the cylindrical system of coordinates (r,θ,z). The mechanicalparameters Young's modulus (F), shear viscosity modulus (η), and density(ρ) are constant in the layer. An analytical solution of the spectralcomponent of the axial displacement in the frequency domain can bederived:

${{Y_{z}\left( {r,z} \right)} = {- {\int\limits_{0}^{\infty}{\alpha^{2}{{J_{0}\left( {\alpha\; r} \right)}\left\lbrack {{A_{1}e^{{- \alpha}\; z}} - {A_{2}e^{\alpha\; z}} + {\alpha\left( {{B_{1}e^{- {\partial z}}} + {B_{2}e^{\partial z}}} \right)}} \right\rbrack}d\;\alpha}}}},{\vartheta = \sqrt{\alpha^{2} - k^{2}}},{k^{2} = {\rho{\overset{\_}{\omega}}^{2}\text{/}\left( {{E\text{/}3} + {i\;\overset{\_}{\omega}\eta}} \right)}},$where J₀ and J₁ are Bessel functions of the order 0 and 1, respectively,and ω is the angular frequency. Four unknown constants A₁, A₂, B₁, andB₂ are defined using four boundary conditions at the layer boundaries.Fixed boundary conditions for displacement on the bottom of the layer;no shear force on the top surface, and a normal force at the focal pointof the transducer are considered. The solution in the time domain wascalculated using the inverse Fourier transform after taking into accountthe duration of the acoustic pulse.

Using the analytical solution of the forward problem, reconstruction ofYoung's modulus and shear viscosity was posed as a minimization problem,i.e. by minimizing the error function defined as the difference betweenthe measured y^(exp) and theoretically calculated displacementsy^(theory) at the point (r=0, z=0):δ=∥y ^(exp) −y ^(theory)(E,η)∥

The density of the lens was assumed to be 1100 kg/m³. To minimize theequation above, a gradient-based iterative procedure was implemented. Inthe minimization procedure, normalized displacement profiles were usedso that only the temporal characteristics of the displacement were takeninto account, not the amplitude of the displacements. This approachavoided the influence of the ultrasound beam attenuation and differencesin acoustic impedance of the materials such as lens and aqueous humor.

Example 5. Uniaxial Mechanical Compression Tests

After the measurements by the US-OCE system, the eye globes werecarefully dissected to extract crystalline lenses for testing with auniaxial mechanical compression testing system (Model 5943, InstronCorp., MA). The lens was centrally positioned between the compressionplates of the device. Prior to the mechanical testing on each lens, a0.004N pre-loading force was applied. The compression speed was set to 2mm/minute. The testing was stopped when the vertical displacementreached 30% the whole thickness. Due to the irregular shape of the lens,it was difficult to directly measure the elasticity based on theconventional compression test method. Thus, an equal-volumetransformation method was adopted to calculate the stress-strainrelationship.

As shown in FIG. 6(a), two parameters of a rabbit lens, the height H,and the maximum diameter D, were recorded before the compression test.The volume of the lens was divided into two regions: the upper region,which is above the x-o-y plane, and the lower region which is below theplane. Each part can be simplified as a half ellipsoid, so the upperregion can be described by the following formula (x²+y²)/(D/2)²+z²/a²=1,while the lower part is described as (x²+y²)/(D/2)²+z²/(H−a)²=1. Thewhole volume of the lens can then be estimated by the integration:

$V = {{V_{1} + V_{2}} = {{{\int_{0}^{a}{{\pi\left( \frac{D}{2} \right)}^{2}\left( {1 - \frac{z^{2}}{a^{2}}} \right){dz}}} + {\int_{0}^{H - a}{{\pi\left( \frac{D}{2} \right)}^{2}\left( {1 - \frac{z^{2}}{\left( {H - a} \right)^{2}}} \right){dz}}}} = {\frac{1}{6}\pi\;{{HD}^{2}.}}}}$

To estimate the stress-strain relationship of the lens, a cylinder withheight H and diameter d, which has the same volume as the lens, wasrequired (FIG. 6(b)). Using the equal-volume equation V=πHd²/4, thecylinder diameter d can be calculated as d=√{square root over (6)}/3DUsing this equation, the stress can be calculated as σ=4F/(πd²)=6F/(πD²)and ε=L/H, respectively. The stress-strain curve was then fitted by theempirical exponent formula:σ=M(e ^(Nε)−1)where M and N are the parameters to be determined. For each deformationcurve, M and N were obtained by using the curve fitting toolbox inMATLAB (Version 2010a, MathWorks Inc., MA). The Young's modulus can becalculated by taking the derivative E=dσ/dε=MNe^(Nε).

Example 6. Results for Lens

The first parameter used to assess the age-related changes inbiomechanical properties of the rabbit lens was the amplitude of thedisplacements as measured by US-OCE. FIG. 7(a) shows the temporaldisplacement profiles of typical young and mature lenses. It is clearfrom either profile that the surface of the crystalline lens starts todeform upon application of the acoustic radiation force. After theremoval of the acoustic radiation force, the surface of the lens beginsto recover to its original position. Based on the high displacementsensitivity of US-OCE, a minimal acoustic radiation force is sufficientto obtain measurable displacements. Here, the displacement is on thescale of micrometers, which helps ensure that the structural andfunctional properties of the crystalline lens are not affected. The datapresented in FIG. 7(b) clearly demonstrate a significant differencebetween the maximum displacements of the young and the mature lenses of3.3±0.1 μm and 1.6±0.4 μm, respectively. It is clear that under the sameexperimental conditions, the maximum displacement of the young lenses isgreater than that of the mature lenses, which is an indicator that themature lenses are stiffer than the young lenses.

The natural frequencies ω=√{square root over (k/m)} in the lenses areshown in FIG. 8. The natural frequency was extracted from the relaxationprocess of the lens surface after the excitation by acoustic radiationforce was removed. The natural frequency is associated with theelasticity of the sample. The recovery process in each of the phaseprofiles was fitted by the model, as shown in FIG. 8(a). During thecurve fitting process, it was found that the damping ratio, ξ, was over0.99, which was approximated to 1 for each eye. Hence, the formulay(t)=A(1+bt)e^(−ωt) was employed to analyze the relaxation process. InFIG. 8(b), the natural frequency, ω, for the two age groups is compared.The natural frequency values of the young and the mature lenses are0.8±0.2 kHz and 2.2±0.5 kHz, respectively. The results show that theyoung lenses had a lower natural frequency than the mature lenses. Inaddition to the maximum displacements, the natural frequenciesdemonstrate that the stiffness of the mature lenses is greater than thatof the young lenses.

The relaxation process is mainly associated with the viscoelasticproperties of the lens. Small oscillations during the recovery processwere observed in both the young and mature lenses, as a result ofdynamic processes in the lens after the rapidly applied force. It shouldbe noted that there is a high variability in the mature samples for boththe elasticity measurements by US-OCE and uniaxial mechanicalcompression testing. This may imply that the effect of age on the lenselasticity varies between individuals.

Quantitative measurements of the mechanical properties of thecrystalline lens based on the US-OCE system required the development ofan appropriate mechanical model and reconstructive procedure. Areconstruction based on the model of a homogeneous viscoelastic layerwas utilized, as presented above. FIG. 9 presents a comparison of thedisplacements measured at the surface of young and mature lenses anddisplacements calculated using the model of a viscoelastic layer withmechanical properties found using the error minimization procedure.Theoretically calculated displacements are calculated for viscoelasticparameters E and η obtained after the minimization procedure: E=2.5 kPaand η=0.37 Pa·s for young lens and E=7.4 kPa and η=0.57 Pa·s for maturelens. Note, that only time characteristics of the displacements wereused and the displacement amplitude was not taken into account. To matchtheoretical and experimental displacement profiles, the acousticpressure at the focal point of the ultrasound transducer in theoreticalmodel was 15 Pa for both for young and mature lenses. While the acousticradiation pressure was not directly measured in experiments, it waslikely the same in all experiments regardless of the animal age.Similarly, the amplitude of the surface pressure was assumed the same intheoretical model. Interestingly, the maximum displacement in the younglens was twice the displacement in the mature lens both in theexperimental and in theoretical results. However, the reconstructionresults show the young lens maximum displacement as three times that ofthe mature lens. In the mature lens, experimental and theoreticalprofiles demonstrate saturation of the displacement at the end of theacoustic pulse, which happens where the elastic response of the mediumcompensated for the acoustic force. These results are in agreement withprevious results, where acoustic radiation force was used to movelaser-induced microbubbles produced in crystalline lenses. However, theUS-OCE method, which combines acoustic radiation force excitation andphase-sensitive OCT measurement, has the advantage of noninvasivedetection. In addition, only a minimal acoustic radiation force isrequired to induce a detectable deformation due to the high displacementsensitivity of phase-resolved OCT detection. This helps to preserve thestructural and functional properties of the ocular tissues.

The result of the reconstruction of Young's modulus and shear viscosityfor young and mature lenses is shown in FIG. 10. While there was asignificant increase in the Young's, the increase of shear viscosity wasless pronounced. These results are in agreement with previous resultsfor bovine lenses, where age-related changes in elastic properties weremore pronounced than changes in viscous properties.

The stress-strain curves and the Young's moduli of the young and themature lenses are compared in FIG. 11. FIG. 11(a) shows typical examplesof the stress-strain curves measured on the young and mature rabbitlenses by uniaxial mechanical compression testing. The stress-straincurves for the young and mature lenses can be fitted asstress=0.1928(e^(12.28×strain)−1) andstress=0.11892(e^(16.83×strain)−1), respectively. At a strain of 0.1,the Young's moduli of the young (n=9) and mature lenses (n=4) are8.2±1.1 kPa and 12.6±1.2 kPa, respectively, as shown in FIG. 11(b). Thisresult clearly shows that the mature lenses are stiffer than the younglenses, confirming the US-OCE results.

As shown in FIG. 9, the results of the model-based calculations are ingood agreement with the experimental data. However, the reconstructedYoung's modulus values are lower than expected based on the results ofuniaxial mechanical compression tests and previous data obtained inbovine lenses. The underestimation in the elastic modulus could be aresult of several factors. The presence of the aqueous humor on thedisplacement of the lens surface was not taken into account, as a freeboundary condition on the surface of the lens was assumed. Thisinfluence may result in the loss of the acoustic energy and Young'smodulus underestimation. In addition, the crystalline lens is aninhomogeneous object, whereas the model used was of a homogeneous layer.Therefore, the reconstruction results may correspond to the cortex ofthe lens, which is softer than the nucleus. Several other factors, suchas the motion of the lens as a whole object, the influence of the lenscapsule, and the spatial distribution of acoustic pressure at the focalzone were not considered in the model, but may also be important.Despite these issues, the current method has been demonstrated as apromising tool for noninvasive assessment of changes in thebiomechanical properties of the crystalline lens in situ. This methodcould successfully applied to other methods of excitation as well asother tissues.

What is claimed is:
 1. A method for measuring biomechanical propertiesof tissues and for differentiating tissue samples using opticalcoherence elastography, comprising: inducing elastic waves in the tissuesamples; detecting properties of the waves using interferometry, lowcoherence interferometry, or optical coherence tomography at measurementpositions along the waves, wherein the measurement positions are thesame or different, and wherein detected properties include measured wavevelocities, measured wave displacement amplitudes, and measured wavedispersion; determining biomechanical properties of the tissue samplesby using the detected properties to create elastic wave displacementprofiles and performing temporal or spatial analysis of the elastic wavedisplacement profiles or one or more of the detected properties of theelastic waves; differentiating the tissue samples having differentelastic wave displacement profiles or one or more different detectedproperties; normalizing the elastic wave displacement profiles for thetissue samples having similar measured wave velocities and needingfurther differentiation to produce normalized wave displacement data;using the normalized wave displacement data to identify tissue sampleshaving faster wave attenuation and slower wave attenuation; andclassifying the tissue samples having faster wave attenuation as tissuesample having increased viscosity and reduced stiffness and the tissuesamples having slower wave attenuation as tissue samples having reducedviscosity and increased stiffness.
 2. The method of claim 1, wherein thetissue samples are ocular or any other soft or hard tissue samples. 3.The method of claim 1, wherein the step of inducing elastic waves is byany internal or external methods such as using an ultrasound/airpuff/laser pulse delivery subsystem.
 4. The method of claim 1, whereinthe step of determining biomechanical properties of the tissue samplesis by calculating Young's modulus using the detected properties of theelastic waves or the elastic wave displacement profiles.
 5. The methodof claim 1, wherein the step of normalizing the elastic wavedisplacement profiles is by dividing the measured wave displacementamplitudes at different measurement positions along the waves by themeasured wave displacement amplitude at an excitation position.
 6. Themethod of claim 1, wherein the step of using the normalized wavedisplacement data to identify tissue samples having faster waveattenuation and slower wave attenuation is by using a customized ratiohaving a formula of:r _((ND) ₁ _(/ND) ₂ ₎=mean(r _(i))±std(r _(i)) wherein$r_{i} = \frac{{ND}_{1i}}{{ND}_{2i}}$ and ND1i and ND2i are normalizeddisplacement data at an i^(th) different measurement position for afirst and a second tissue sample, wherein if r_((ND1/ND2)) issignificantly greater than 1, the second tissue sample is identified ashaving faster wave attenuation and the first tissue sample is identifiedas having slower wave attenuation, and wherein if r is significantlyless than 1, the first tissue sample is identified as having faster waveattenuation and the second tissue sample is identified as having slowerwave attenuation.
 7. The method of claim 6, wherein the step of usingthe normalized wave displacement data to identify tissue samples havingfaster wave attenuation and slower wave attenuation is repeated fordifferent tissue samples.